Non-Markovian Models of Blocking in Concurrent and Counter cur rent Flows 
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We investigate models in which blocking can interrupt a particulate flow process at any time. 
Filtration, and flow in micro/nano-channels and traffic flow are examples of such processes. We 
first consider concurrent flow models where particles enter a channel randomly. If at any time 
two particles are simultaneously present in the channel, failure occurs. The key quantities are the 
survival probability and the distribution of the number of particles that pass before failure. We then 
consider a counterflow model with two opposing Poisson streams. There is no restriction on the 
number of particles passing in the same direction, but blockage occurs if, at any time, two opposing 
particles are simultaneously present in the passage. 



Introduction. Processes involving the flow of particles 
through channels may entail blocking or failure. A good 
example for a concurrent flow is provided by the indus- 
trially important process of filtration [lH5]. In particular, 
the model of Roussel et al. [S] successfully accounted for 
experimental data by assuming that clogging may occur 
when two grains are simultaneously present in the vicin- 
ity of a mesh hole, even though isolated grains are small 
enough to pass through the holes. A conceptually simi- 
lar situation is a flimsy bridge that can only support the 
weight of one car at a time. If ever two cars are on the 
bridge at the same time, it collapses. 

A second class of processes involves two counterflowing 
streams of particles. For example, in remote areas many 
of the roads are single-track. Two approaching vehicles 
cannot pass each other except at rather infrequent, and 
short, passing places. In this situation we would like to 
know the failure probability of finding two opposing cars 
in the stretch of road between two passing places. 

Many traffic models based on lattice gases have been 
proposed [6lfl3] including the totally asymmetric simple 
exclusion processes (TASEP) [TH [TS] and related mod- 
els ITS'. The so-called bridge models [T71 - f2^ consider 
two TASEP processes with oppositely directed flows, but 
allow exchange of particles on the bridge. 

Similar processes are also found in numerous biological 
applications involving channels. Examples include bidi- 
rectional macromolecular flow in microchannels [23], ion 
channels that can be clogged by toxins or medicines [24l - 
126] ■ and the antibiotic gramicidin that forms univalent 
cation-selective channels of 0.4nm diameter in phospho- 
lipid bilayer membranes. The transport of ions and water 
throughout most of the channel length is by a single file 
process; that is, cations and water molecules cannot pass 
each other within the channel [27] . 

In this Letter we propose, and obtain exact solutions 
for, stochastic models in which particle fiow in a channel 
can be instantaneously interrupted by a clogging event. 
The quantities of interest are the probability of block- 



age (failure) as a function of time and the final outcome, 
i.e. the number and type of particles that get through 
the channel before blockage occurs. These models are 
complementary to the lattice gas models in that they are 
continuous in both space and time and are most appro- 
priate for low density flows. 
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FIG. 1. (a) Concurrent flow model: Particles enter the left 
hand side of a channel of length L at a (mean) rate A. Block- 
age occurs when two particles are simultaneously present in 
the channel, (b) Counterflow model: Two opposing streams 
of particles enter the left and right hand sides of the channel at 
rates Ai and A2, respectively. Blockage occurs whenever two 
opposing particles are simultaneously present in the channel 
but there is no constraint on the number of particles moving 
in the same direction. 

Concurrent flow model. Particles enter a passage of 
length L according to a homogeneous Poisson process 
where 

^„W = ^exp(-At) (1) 

gives the probability that n particles enter the passage 
in the time interval (0,t). We assume that all particles 
move with constant velocity v so that the transit time, 
T — L/v, is constant. Blockage (failure) occurs at the 
instant when two particles are present in the channel at 
the same time (see Fig.[l]a). This leads us to consider the 
survival probability Ps{t), the probability that blockage 
(failure) does not occur in the time interval (0, t). Clearly, 
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Ps(0) = 1 and Ps{oo) = 0. The probability that blocking 
occurs between time t and t + dt is given by f{t)dt where 

fit) = -dpsit)/dt. 

To solve the model we introduce the n particle survival 
probability qs{n,t) which denotes the joint probability 
of surviving up to t and that n particles have entered 
the passage during this time. The survival probability 
is simply Ps(t) = X]n>o evolution of the 
qsin,t) is given by: 

^^i^=Xqs{n-l,t-T)e-^^-Xqs{n,t), n>2 

(2) 

The final equation implements a non-Markovian con- 
straint: in passing from the state {"not blocked", n — 1} 
to the state {"not blocked", n} it is necessary that no 
particle enter in the previous time interval {t — r, t) and 
that a single particle enter in the time interval t to t + dt. 
These probabilities are given by e"'*''^ and Xdt, respec- 
tively. By introducing the generating function we obtain, 
as detailed in the Supplementary Material (SM), 



{X{t-nT)l_ 
' {n + l)\ 

where 9{t) is the Heaviside function. 
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-At 



(3) 
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FIG. 2. Probability that blockage occurs as a function of time 
for Xt — 2, 1,0.5 (most peaked to least peaked). Cusps are 
present at t/r — 1 (vertical dotted line). The inset shows the 
survival probability for the same parameters together with 
the asymptotic approximation, Eq. Q (dashed lines). 

The long time behavior of the survival probability 
can be obtained by approximating the sum in Eq. ^ 
as an integral and evaluating it using the Saddle-Point 
(Laplace) method. The result is 



Ps{t) 



(4) 



where Lw(x) is the Lambert-W function (see Sect. 1.2 of 
the SM). From its small x behavior one can deduce that, 
when Ar ^ 1, the exponent of the exponential decay. 
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FIG. 3. Distribution of the number of particles that pass be- 
fore blockage occurs for different values of Xt (left and middle 
figures) . The right graph shows the average number as a func- 
tion of the same parameter. 



A — Lw{\t) ^ X^T, depends nonlinearly on the rate A. 
This complexity arises from the large possible number of 
event sequences before failure. 

The mean survival time is given by 



{t) = / Ps{t)dt = 



2e' 



1 



A(e^^ - 1) 



(5) 



and is consistent with Eq. (|4| when At ^ 1. 

Figure [2] illustrates the time dependent properties of 
the concurrent flow model. The curves showing the prob- 
ability of failure at time t, f{t), exhibit a cusp at t = t. 
This is the first time at which particles that have en- 
tered previously can exit the channel (which is empty 
at i = 0), leading to a rapid decrease in the probabil- 
ity of blockage. The intensity of the cusp depends on A, 
{df /dtl^.- — df/dt\T.+ = X^e"^^), and is less pronounced 
for large A as a second particle is more likely to enter 
soon after the first, causing blockage. The inset shows 
the survival probability and confirms the accuracy of the 
asymptotic expression, Eq. (|4|. 

A further quantity of interest is the distribution of 
number of particles that exit the channel before blockage 
occurs. If n particles have entered the passage at failure, 
the number that have successfully traversed is m = n — 2. 
At least two particles must enter before failure can occur. 
Let h{m) denote the probability that when failure occurs 
m particles have exited. If At^ denotes the time interval 
between the entry of the i*'' and the {i + 1)*'* particle, 
then 



h(m) 



[]Pr(At, 



Pr(AWi < r) (6) 



Using that Pr(Ai > r) = e"-^^ and Pr(At < r) = 1 - 



we find 



/i(m) = e~'"^^(l-e-^^) 



(7) 



The most probable situation is that no particles pass 
before failure for all values of A. The mean number that 
pass before failure is 



1 



e^^ - f 



(8) 
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which has the expected asymptotic behavior: (m) — 
e~^'^ for Ar large and (m) — >■ (Ar)"-'^ for At small. Figure 
[3] illustrates that with decreasing At the difference be- 
tween the mean, (m) , of m and its most probable value 
(always 0) increases and h{m) becomes flatter. 

The above results can be generalized for a distribution 
of transit times. Hipij) is the normalized distribution of 
transit times and assuming that A is constant, the mean 
survival time is 



{t) = 



2-^(A) 



(9) 



A(1-^(A)) 

and the mean number of particles that pass before failure 



IS 



1 - ^(A) 



(10) 



where tilde denotes the Laplace transform, 



r 
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-At 



1p{T)dT. 



Counterflow model. In this model particles of type 1 
enter the channel of length L at the left at a rate Ai and 
move towards the right at a speed vi. Particles of type 
2 enter at the right at a rate A2 and move to the left at 
speed V2'- (see Fig. [T]3) The transit times arc ti = L/vi 
and T2 = Llv2, respectively. We assume that particles 
enter according to a Poissonian distribution so that the 
probability that ni (712) particles of type 1 (2) enter the 
left (right) side in the interval (0, t) is given by: 

(A.t)"' 



-e~^'* within 1,2 



(11) 



A blockage occurs if, at any time, particles of both species 
are present in the channel. Before this situation arises an 
arbitrary number of particles can transit the passage in 
both directions. If A^t; < 1 the average time interval 
between entry of particles of type i is longer than the 
transit time. If, on the other hand, A^t^ > 1, a backlog 
of particles of type i is likely to be present. Thus the 
former situation is more relevant physically. 

We now outline the solution method. The device that 
allows us to obtain an analytical solution in this case is 
the introduction of functions pk{ni, 71,2; t) that denote the 
probability that the system has survived until time t and 
nil particles of type 1 and 712 of type 2 have entered the 
passage and the last particle to enter the passage was of 
type k = 1,2. This choice provides a complete parti- 
tion of the event space into disjoint events allowing us to 
write p{ni,n2;t) = Pi{ni,n2\t)+p2{ni,n2]t),ni,n2 ^ 0, 
Ps(0,0;f) = pi(0,0;i) = p2(0,0;t) (by convention) and 
vX^) ^iln,=QY.n^=oV{ni,n2;t). ^ 

The equations describing the time evolution of the 
probabilities Pk{ni, n2', t) are e.g. 

dpi{ni,n2\t) 



dt 



-(Ai + \2)pi{ni,n2]t) 



+ Ai[pi(ni - I,n2;t)+P2(ni - 1, n2; t - T2)e-(^i+^^)"^] 

(12) 



for ni > and n2 > 0. The last term of this equation im- 
plements the constraint for the " not-blocked" state of the 
channel (see Supplementary material): in passing from 
the state {"not-blocked", ni — 1, 712, last particle entered 
— type 2} to {"not-blocked", ni,n2, last particle entered 
= type 1} it is necessary that in the previous time interval 
(t — T2, t) (i) no particle of type 1 or 2 enters the channel 
(given by the exponential term) and (ii) a single particle 
of type 1 enters between t and t + dt. In analogy with 
Eq. ([2]), this is indicative of the non-Markovian nature 
of the process. The evolution equation for P2{ni, n2] t) is 



obtained from Eq. (12 1 by symmetry. 

In addition, we have to consider the time evolu- 
tion of the "boundaries" (0,712) and (ni,0): obviously 
P2("i,0;i) =pi(0,7i2;t) = for 77,1,712 > 1. 



dpi{ni,Q;t) 
dt 



= -(Ai + A2)pi(77i, 0; t) + \ipi{ni - 1, 0; t) 

(13) 

with Til > 1 and a corresponding equation for P2(0, n2]t). 
To complete the configuration space, one must introduce 
the probability that no particle is created in the time 
interval (0,t), p(0, 0;t) and one has^^^^^^ = — (Ai + 
A2)p(0,0;t) with p(0,0;0) = 1 with solution p(0,0;t) = 

g-(Ai + A2)t_ 

As for the previous model, the solution is obtained by 
introducing a generating function: 

00 00 

G(zi,Z2;t) ^ ^ 22X771,712;*) (14) 
m— 712—0 

from which the survival probabiHty can be found as 

p,(i) = G(l,l;0. 

After some calculation (see Sect. 2.1 of the SM), we 
obtain 



G(zi,Z2,7i) = 



1 



1 + Ai + A2 
A2Z2 



Xizi 



Ti + A2 + Ai(l — Zi) 

A1A2Z1Z2 
A 



M + Ai+A2(l-Z2) 

|^g-(Al+A2+il)Tl _|_ g-(Ai-|-A2-|-u)r2 _|_ g-(Ai+A2+u)(Tl+T2) 

Ai^i ^ A2Z2 

77-|-A2-f'Ai(l-Zi) ?7 4-Ai+A2(l-Z2) 

(15) 

where A (u + A2 + Ai(l - zi)){u + \i + A2(l - Z2)) - 
AiA22:iZ2e"^^i+^=+"^^^i+^"^ This is the principal result 
for the counterflow model from which most properties of 
interest can now be easily obtained. In particular, the 
mean survival time, {t) = ps{u ^ 0), is 



t) = 



1 



+ 



Ai + A2 
Af + Al 
A1A2 



1 



p-(Al+A2)Ti _j_ g-(Ai+A2)T 

I - e-(Ai + A2)(Ti+r2) 

1 



,-(Ai+A2)(ri+T2) 



(16) 
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The three contributions have a simple physical mean- 
ing: the first term corresponds to the situation where 
no species exit the passage before failure. The second 
term corresponds to situations where an even number of 
changes of species occurs before failure, the last term to 
the situations where an odd number of changes (larger 
than 1) of species occurs before failure. 
Note that when (Ai + X2)Ti « 1, 



1 



AiA2(ri + T2) 



(17) 



corresponding to a regime where a large number of event 
sequences contribute to the survival probability. 




FIG. 4. Probability of blockage, f{t) — —dps{t)/dt, as a 
function of time for the counterflow model. Ai = 0.5, A2 = 
0.2,T2/ri = 1 (red), T2/tx = 2 (blue), T2/t± = 3 (green), 
r2/ri = 5 (black) 




FIG. 5. Probability of different outcomes in the counterflow 
model as a function of pi = Ai/(Ai + A2), Ai + A2 = 1 and 
Ti = r2 = 1. The solid curves indicate the probability that 
zero (= /i(0,0)), one (= /i(l, 0) + /i(0, 1)) and two (= h{2,Q) + 
h(0,2) + h{l, 1)) (top to bottom) particles exit the channel 
before blockage occurs. 



It is more difficult to obtain /i(mi, TO2), the probability 
that mi particles of type 1 and m2 particles of type 2 
exit the passage before blockage occurs. This is because 
there is no simple relationship between mi,m2 and the 
numbers ni,n2 that have entered the passage as is the 
case in the concurrent flow model. However, the first few 
may be obtained by direct calculation: 



/i(0,0) =pi(l- 



— A2ri 



(20) 



It is possible to perform a term-by-term inversion of the 
Laplace transform to obtain the time dependent survival 
probability. For the case Ai = A2 = A, ri = r2 = r the 
result is (see the SM) 



Psit) 



-2\t 



k=0 



_g-2At g-A(t+2fcT) 



2k 



5^^^(A(t-2fcr))' 



i=0 



2j2e{t- {2k + l)T) 



k=0 



i-iy 



2k+l 

^A(t+(2fc+l)r) 

1=0 



(A(t-(2fc + l)T)' 
(18)' 



For a given time, the solution contains a finite number 
of nonzero terms. At large t, by using Laplace's method, 
we obtain 



Ps{t) 



(19) 



Note that when At ^ 1, ps(i) ~ e "^^ consistent with 
Eq. (17). The average survival time is dominated by this 



regime. The general solution is given in the SM and Fig. 
|4] illustrates a particular case. Note the presence of two 
cusps (ti 7^ T2) corresponding to the two transit times. 



/l(mi,0) =^7*1^2(6"^"^^ _ g-(Ai+A2)rig-Air2>,^ > ^ 

(21) 



— n T)"2 „ p-(Al+A2)T2 -A2T1 



/i(0,m2) =PiP2^{e 



), ma > 1 
(22) 



and 



/i(l, 1) = piP2e 
+ P2{e^'^' 



(Al+A2)(Ti+T2)|-p^|-gA2T2 
-A1T2 \1 



-A2T1 \ 



— e 



(23) 



where Pi — Ai/(Ai -I- A2), « = 1,2. See the SM for details. 
The behavior of these functions is illustrated in Fig. [5] 
for the non-restrictive situation of a constant total flux 
Ai + A2 = 1 apportioned continuously between the left 
and right hand streams. As in the concurrent flow model, 
the most likely result is that blockage occurs before any 
particles exit. 

We conclude with an illustration of the theory: Let us 
suppose that a single-track road is O.bkm long and on 
average 10 cars enter each side per hour. If we further 
assume that all cars travel at a constant speed of 50km /h 
then the survival probability after 5 minutes is 0.876 and 
after 30 minutes it is 0.436. The mean survival time is 
36 minutes. 

In summary we have developed stochastic models to 
describe the probability of blocking in diverse physical 
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applications involving particulate flow. Both models can 
serve as the starting point for more refined models tai- 
lored to specific applications. For a filter composed of M 
independent channels, the fraction that is active at time 
t is just Mps{t). With more effort, connected channels 
and reversible blocking can also be treated within the 
same framework. Clustering of the particulate streams 
can be modeled using an inhomogeneous Poisson process 
where the intensity is time-dependent, A(t). 

We thank O. Benichou, G. Oshanin, Ch. Pouzat, L. 
Signon and G. Tarjus for useful discussions. 
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